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Abstract 

We present a two-dimensional spring-block model of earthquakes with the full nonlinear 
equations for the forces in terms of displacements. Correspondence of our linearized version 
to earlier models reveals in them inherently asymmetric elastic properties. Our model gener- 
alizes those models and allows an investigation of the effects of internal strains and vectorial 
forces. The former is found to be relevant to critical properties such as that described by 
the Gutenberg-Richter law, but the latter as well as nonlinearities up to second order in 
displacements are irrelevant. 
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In 1967, Burridge and Knopoff[|l] introduced a one-dimensional (ID) system of springs 
and blocks to study the role of friction along a fault in the propagation of an earthquake. 
Since then, many other researchers investigated similar dynamical models of many-body 
systems with friction, ranging from propagation and rupture in earthquakes [2-10] to the 



fracture of overlayers on a rough substrate|ll 



Among these developments, a deterministic version of the ID Burridge-Knopoff (BK) 
model was analyzed by Carlson and Langerp| and the same model but in a quasi-static limit 
was studied by Nakanishi0. A 2D quasi-static variant was first simulated by Otsuka[Q and 
later by Brown, Scholz and Rundle^, who introduced a discrete version that was formulated 
as a cellular automaton. A similar model with continuous local variables, generalizing the 
model of Bak, Tang and Wiesenfeld^, was studied by Olami, Feder and Christensen[]7| 
(OFC). Contrary to previous models, the model is non-conservative and derivable from a 
2D BK model (under certain limit, see below). It produces robust critical behavior which 
depends on the conservation level. OFC argued that such dependence explains the variance 
of the exponent in the Gutenberg- Richter [p!2| law observed in real earthquakes |T3| . 



Model. — As before, our model consists of a 2D array of interconnected blocks in contact 
with a rough surface. Each block is also connected to a rigid driving plate by a leaf spring 
whose spring constant is Kl. The plate moves at constant, infinitesimal speed. The coor- 
dinates for the attachment of the leaf springs on the moving plate, labeled by where 
^ ^ h j ^ L, form a square lattice with lattice constants ai, 02 (hereafter the subscripts '1' 
and '2' refer to properties in the x- and y-direction, respectively). Driving the plate thereby 
induces stress between the array and the plate. Each block in the bulk is connected to four 
nearest neighbors via coil springs whose spring constants are Ki and K2 and unstretched 
lengths h and l2- We restrict ourselves to the situation where a > I, and the displacements 
Xij, Uij measured from (z,j) fulfill Xjj <^ ai and <^ 02, so that Hooke's law applies. 
When the net force on a block exceeds a static frictional threshold Fth-, the block slips 
instantaneously to a new equilibrium position. 
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While this is an obvious extension of the ID BK model[||] to two dimensions, equations 
from previous models ^, fail to describe it correctly. The reason is that elements in 
those models have been restricted to move in one direction only. To our knowledge, this 
restriction has not been justified beyond the reason of simplicity. To address its relevance, 
we lift this limitation and start with the proper nonlinear equations for forces as functions 
of the displacements of blocks. The net force Fi^ (which equals the friction due to the rough 
surface) on the block at (i, j) is a vector: 

Fi,j = fL + f{i+l,j)-{i,j) + f{i-l,j)-ii,j) + + (1) 

where f(i±i,j±i)-{ij) is the force exerted by a neighboring block and fi the loading force by 
the driving plate. For example, we have 

J{i+i,j)-(i,j) = J^iiai + Xi+ij - Xij ==j, [Z) 

= Mvw-Vi,- I =]• (3) 

[ai + Xi+ij - XijY + {yi+i,j - VijY 



Then the x component of the force in Eq. (|l|) takes the form: 

F^. = fl-Ki [2xij - Xi+ij - Xi-ij + 

(ai + Xi^ij — Xij)li (ai — Xi-ij + Xij)li 



(ai + Xj+ij - XijY + (z/i+i.i - VidY \J l*^! ~ ^i-iJ + ^ijY + iVij ~ yi-idf 



(4) 



\/ (^2 + yi,j+i - IJijy + - XijY \J («2 - Vij-i + Z/ij)^ + - a^ij-i) 

and, by symmetry, Ff^ follows by switching x ^ y, i j , li ^ I2, cli ^ a2 and Ki ^ K2. 

Linear version. — In order to compare with the OFC model, we expand Fi j to first 
order in x's and y's. Specifying fi = {—KiXij,0) as in 0, one can easily show that a 
slip of the block at results in the following force redistribution |jT5|, p!6[| : 



Fr,±i - F:^,±i + S2aa,Fl^, (5) 



^i±i,j - 




a 


^i,j±i - 


^ T?y 




p. . _ 


0, 





where = (ai — li)/ai is the internal strain in the x direction, and similarly for 5*2. 
cr = K2/K1 and K = Kl/Ki are measures of anisotropies in the couplings. In the bulk, ai 
and 0^2 are given by: 

" 2(1 + S^a) + ' 2(l + 5i/a)- 

Their values at boundary sites depend on the boundary conditions (see below). 

Using Eq. (5), our model can be described as a coupled map lattice as was done in 0. 
If we impose = at all sites as in H, 0, so that (to first order) all Ffj = 0, we readily 
recover the OFC model in the limit S2 —* 1, i.e., for maximal internal strain along the y 
direction. This also follows from the general relation Eq. @j, as becomes linear in x for 
5*2 = 1 and y = 0. The physical rationale behind this correspondence is that in [^, |^, 0, 
both the loading springs and those in the array along y were chosen as leaf springs, which 
by definition cannot be extended along their length but have a restoring force linear in the 
transverse displacement. They act like fully stretched coil springs. With leaf springs, the 
forces become scalar and the equations linear. However, it introduces an artificial asymmetry 
into the system, which to our knowledge has no analog in nature. This manifests itself most 
clearly in the elastic moduli P^, For the OFC model and the like, Cim = Kiai/a2, 
C2222 = 00, C1122 = 0, and Ci2i2(shearmodulus) = K2a2/ai. Clearly, setting Ki = K2 does 
not render the model symmetric. For our model, the asymmetry is absent: Cuu = Kiai/a2, 
C2222 = ^^202/01, C1122 = 0, and C1212 = {KiK2SiS2aia2) / {KiSial + K2S2al). 

For general strain (0 < 5*2 < 1) and = 0, the linear version of our model coincides with 
the anisotropic OFC model0 with our 5*20- corresponding to OFC's K2/K1. Consequently, 
except in the unphysical limit 5*2 ^ 1, asymmetric loading on a symmetric model {k = 
a = 1, S2 = Si) gives rise to intrinsically asymmetric force redistributions. The statement 
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regarding the variance in the Gutenberg- Richter law as a result of different elastic parameters 
has to be reinterpreted, in the present context, as a result of variances in both internal strains 
and elastic parameters. 

Allowing for i/ij ^ 0, Eqs.(5,|D yield the total change of force after a block in the bulk at 
slips by a distance {6x,6y): 

SF"" = -Kl6x = -KaiFr'j ; = 0. (7) 

The change on the boundary depends on the boundary conditions, which are defined by 
the number of neighbors a boundary block is attached, and the way it is loaded. Following 
OFC, in "free" boundary conditions (FBC) the sites on edges have three neighbors while 
those at corners have two. In "open" boundary conditions (OBC) all sites has four neighbors. 
Both have uniform loading throughout the system by an amount AF^ = KiAx, where Ax 
denotes the displacement of the driving plate relative to the array. 

For FBC, the a's on the boundary differ from the bulk values in Eq. (P): 

{x=l,L) 1 (y=l,L) 1 (corner) 



(x=l,L) _ 1 (y=l,L) _ 1 (corner) _ 1 \.°) 

"2 2+51/(7' "2 1+251/(7' "2 l+5i/cr- 

It turns out then that Eq. (|^) holds at all sites. While the change of F^ is balanced on the 
average by loading, the spatial mean F^ is always zero due to Newton's third law and the 



fact that the system is isolated for FBC. Since a slip is decided by y (F^)^ + (F^)^ — Ffh, the 
two components are locally coupled in the dynamics. We find in the steady state that the 



fluctuation of F^ tends to zero, so that this coupling to a conservative field is irrelevant |T7 
i.e., the OFC limit is recovered. 

For OBC, Eq. (D holds at all sites but Eq. (0) is modified to SFy/F^ < at boundary 
sites due to their coupling to the external imaginary blocks. This implies that a system with 
arbitrary initial spatial mean, F^, will always flow in steady state to the F^ = (stable) 
fixed point (i.e., the OFC model). We conclude, for OBC too, that the steady state of our 
linearized version is the same as that of the OFC model. 
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We now remark on the boundary conditions. We first reveal a flaw in previous treat- 
ments of OBC[0. Effectively, the ' 'imaginary layer of blocks" in OBC corresponds to a rigid 
frame that attaches to the array by springs. Spatially uniform loading in Ref. implies 
no relative displacement between the frame and the array during loading. Since the imag- 
inary layer never slips, it is permanently pinned on the rough surface. Although this is an 
unphysical setup, it has not been appreciated before in ||^, Q because it is not obvious in 
the Fij conflgurations. Working with the linearly related Xij instead, we show in Fig. 1 
how the initially square array of blocks is distorted due to the pinning. At long times, 
the distortion can be arbitrarily large so that there is no meaningful steady state for this 
model. To make physical sense, the frame has to move along with the driving plate. It is 
then intuitively clear that loading cannot be uniform: due to pulling and pushing by the 
frame, boundary blocks are loaded more than in the bulk (e.g., AFj^^ = (1 + k^^) /S.F^^y^ . 
Incorporating such non-uniformity, we observe numerically |ll5[| that the model becomes non- 
critical and reaches periodic states much like with periodic boundary conditions0. Second, 
for FBC, we implicitly assume in Eqs. (5) and that the internal strains are somehow 
sustained. Otherwise, due to missing neighbors at the boundary, the equilibrium positions 
of blocks would be shifted inward relative to a regular lattice. In practice, the locations of 
loading springs on the driving plate may be adjusted to eliminate the shifts. In any case, 
bulk properties are not affected because the shifts decay exponentially into the bulk over a 
distance[=l/ ln(2 + K)] of order unity |I5|. The use of leaf springs along y direction in [^, ||, |^ 
may also be interpreted as another means of imposing the strains. 

Nonlinear version. — Eq. (^) allows for an investigation of nonlinear effects which 
naturally arise from the full forcing. However, it is no longer possible to describe the model 
as a coupled map lattice [cf. Eq. (5)]. Instead, one has to keep track of the displacements in 
order to determine whether \F\ > Fth — if so, one has to solve for the displacements (x, y) 
that deflnes the zero-force position of the block in terms of its environment. To n^^ order. 
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this means to solve two coupled equations of the form [for brievity, q = 

"19 d 
Fq = = fl - E ^^p,Mp,i - E K<? + Vp,,TrrGp,,{x, 2/)|(.,j,)=(o,o)}, (9) 

p m=o ^-^ 

plus a corresponding equation [F^ = 0) obtained by the transformation after Eq. (^. Here 
p denotes the four nearest neighbors (z ± 1, j ± 1), and the symbols stand for: 

Kp,, = K2, Gp,,ix,y)= \ for p = (z, j ± 1); (10) 

Up,q = Xq-Xp, Vp,q = yq-yp, fOT p = {t ± 1 , j) , {i , J ± 1) . 

Note that we have implicitly assumed that no nearest neighbors can be supercritical at the 
same time, due to infinitesimal loading rate. 

Apparently this program is not computationally efficient for large n. But n = 2 is simple: 
It is easy to see that all second order terms {x'^,y'^,xy) cancel in the bulk (i.e., when 
has four nearest neighbors), resulting in two coupled first order equations for {x, y) in terms 
of F and (x, y) before slipping: 

= {xij - Xij)[a^^ + -^{yij+i - Vij-i)] 

+iyi,j - yi,j)[^iyt+i,j - Vi-ij) + ^{xij+i - Xij.i)], (ll) 

plus a corresponding equation for Ffy Note that it yields nonlinear dependence of (x,y) on 
(x, ?/) via those terms proportional to l/a^. 

Based on Eq. ([TT|) , we have performed simulations for different boundary conditions. Due 
to missing neighbors of the blocks on the edges for FBC, second order terms in (x, y) survive 
and one has to solve higher order (3''*^ and 4**^) equations for the equilibrium positions of the 
boundary blocks. Remarkably, even with substantial nonlinearities (determined by l/o?)-, in 
none of the cases of FBC, OBC and PBC do we find the critical behavior to deviate from 
that of the linear cases. 

To summarize, we have investigated the critical behavior of a general 2D spring-block 
model of earthquakes, within the context of the BK model under the quasi-static limit. We 



discover that internal strains are important ingredients in the variance of the Gutenberg- 
Richter law. Nonlinear perturbations, to which the OFC model displays remarkable stability, 
are not introduced arbitrarily but generated generically at the microscopic level from the 
elastic spring-block interactions. Among possible extensions, the case of internal compres- 
sional strains (i.e., S < 0), and of different loadings, including those along both directions 
[i.e., fi = —{Kf^x, Kf^y)] and spatially non- uniform ones (e.g., shear applied through the 
boundary), are being pursued[|l^]. 
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Figure captions. 



Fig. 1. Snapshot of a configuration of blocks after 2000 avalanches, for L = 20, ai = 02 = 1, 
^1 = ^2 = 0.1 and a = K = Ki = Fth = 1, showing the unphysical effect of pinned frame 
(denoted by □) in OBC. The array is pulled to the right. The system evolves according to 
Eq. (0). 
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